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Abstract 


A brief review of Stokes' problem for the ellipsoid as 
a reference surface is given. Another solution of the problem 
using an ellipsoidal kernel, which represents an Iterative 
form of Stokes' integral, is suggested with a relative error 
of the order of the flattening. On studying of Rapp's method 
in detail the procedures of improving its convergence are 
discussed. 
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Introduction 


The determination of the shape of the global earth Is 
one of the topics pursued in physical geodesy, while the Inves- 
tigation of the earth's shape consists essentially in eval- 
uating the height anomaly and the deflection of the vertical 
by using the gravity anomaly Information at every point on 
the surface of the earth. Since 1849 Stokes developed a formula, 
which states the disturbing potential of the earth can be 
derived from gravity ancwaaly values at the geoid, his famous 
formula has been traditionally regarded as a basis to determine 
the earth's figure and extensively applied to practice. But 
there are two important limitations of its practical appli- 
cation in geodesy. That is to say this formula is valid only 
for (a) spherical boundary surface, and (b) distribution 
of gravity values over the geoid instead of the topographical 
surface. With spherical approximation it implied that the 
error of the order of f , the flattening of the earth, will be 
certainly introduced into the computed results. For instance, 
the absolute effect of this relative error on the geoid undulation 
N reached fN s 0.7 meters when the geoid undulation is 100 
meters (Fang Tsiin, 1973, p.236). Although this error is usually 
permissible for normal purposes, it should be taken into consid- 
eration for special investigation like the accurate determination 
of the figure and gravity field of the earth. In addition, 
the requirement that the gravity values must refers to the 
geoid and no masses can lie outside the geoid is so stern 
that it is impossible for this requirement to be exactly satisfied 
even by using various calculation procedures to reduce gravity 
and regularize the surface of the earth's geoid. Therefore, 
from both theoretical and practical points of view, we cannot 
use the classical equation of Stokes for analysing in detail 
the earth's geoid; instead wc have to develop modern theories 
and techniques. During the last 40 years or so there has 
been a great amount of work being done in this field aiming 
to improve and develop the classical Stokes' theory. It should 
be pointed out that the most important contributions have 
been made by Molodensky (1962'', Zagrebin (1956, 1965), and 
Pellinen (1981). They approached the problem in two different 
ways. Molodensky gave up the traditional concept of the geoid 
of the earth and, instead. Introduced an auxiliary geometric 
surface, the so-called telluroid, to reduce the problem to 
the solution of Laplace boundary value problem of mixed type 
on the telluroid, which is called the geodetic boundary value 
problem in geodesy. Thus, we can find the height anomaly 
and the deflection of the vertical directly by using gravity 
anomaly data measured on the actual surface of the earth. 

Zagrebin sought a solution of the boundary value problem with 
the ellipsoid as a reference surface to carry out the investi- 
gation of the gravity field of the regularized earth. Both 
these theories establish the most elaborate form of the so- 
called model approach to physical geodesy. 


- 1 - 


I • • 


In the present report, we are not going to discuss 
Molodensky's problem, which can be found in the authoritative 
work ADVANCED PHYSICAL GEODESY (Moritz, 1980); Instead we 
are just planning to try to given a brief review of the state 
of art about the Zagrebln's problem, the problem of Stokes 
for the ellipsoid. Further-more, we will place stress on 
th' new techniques which are appropriate for practical appli- 
cation, such as those developed by Rapp (1981). In addition, 
we also contribute another solution with iterative form for 
the problem of Zagrebin using an ellipsoidal kernel, which is 
similar to Molodensky's (1962, p.53) and Koch's (1968, p.22). 
Moreover, on studying Rapp's method in detail, we put 
forward some new designs for improving the convergence of 
his method. 


2. The Form of the Solution 


Since the Zagrebln's solution was published a large number 
of authors have been devoted to solving the Stokes' problem 
for a reference ellipsoid with a relative error of the order 
of the flattening. Looking into all works it should be found 
that all available approaches, from the mathematical point 
of view, could be classified into two groups: the direct 

solution - series approach, expanding the anomalous potential 
of the earth into a series of harmonics, such as Zagrebln's 
(1956, 1965) and Bjerhammar's (1962), and the indirect solution - 
ellipsoidal kernel, reducing this problem to the solution 
of the Fredholm's integral equation of the 2nd kind by means 
of either the third Green's formula or the theory of the poten- 
tial of surface layer, such as Molodensky's (1962), Fang Tsun's 
(1965), Koch's (1968) and Lelgemann's (1970), etc. 

In spite of using any methods in these two groups to 
solve the problem of Stokes for an ellipsoid of revolution 
the solution finally obtained could be written as the fol- 
lowing form: 

(1) Tg = T + 6T 

or from the Bruns' formula, 

(2) Ng = N + 6N 

where T or N is the principal part of the solution, expressed 
by Stokes' equation, which characterizes the fundamental pro- 
perties of the gravity field of the earth, and 6T or 6N 
is a correction terra, which responds to the compression of 
the earth. 

In case of the direct approaches, the solutions are rep- 
resented as an open form, the series of harmonics; otherwise, 
a closed form, the Stokes' integral. 


- 2 - 




We can see that the series treatments are not only tedious 
in derivation, but also very complicated in final form obtained. 
For example, for the Zagrebin's solution it is necessary to 
carry out not less than nine Integrals over the reference 
surface in order to obtain Ng of one point. Thus it can 
be seen that such a kind of solution is not appropriate for 
practical application. In view of this situation, another 
way, using an ellipsoidal kernel, to approach the problem 
was suggested by Molodensky, and afterwards developed by Fang 
Tsun, and Koch, etc. From the pure mathematical point of 
view, the solution of ellipsoidal kernel methods, which is 
expressed as a set of iterative Integral forms, is graceful 
and strict, •'nd convenient for theoretical analysis. However, 
unfortunately, it should be troublesome to employ them to 
evaluate geoid undulations through potential coefficient infor- 
msition. For this reason, it is necessary to Improve and develop 
the kernel methods so as to suit them to the need of the modern 
practice in geodesy. Before discussing this point, we are 
going to contribute another solution of the Zagrebin's problem 
with the aid of ellipsoidal kernel in the following section. 


3. Another Solution of the Problem Using an Ellipsoidal Kernel 


Most authors utilize the third Green's theorem to seek 
the solution of the Zagrebin's problem; Instead we still 
start with the theory of the potential of a simple layer as 
Koch (1968) did, but a special kind of coordinate transformation 
io employed, by which the whole procedures of derivation is 
significantly simplified. Because of simplification , the approach 
developed here should be meaningful. Now, we will show the 
derivation as follows. 


Assume the eqauation of the reference surface is given 
by 


(3) 




1 


First of all, let us introduce the coordinate transfor- 
mation 


A 


(x,y,z) 




(4) X =» E coshy cos 6 cos 1 
y * E coshy cos 6 sin X 
z » E sinhy sin 6 

E ® (a^ - c^ “ c e 
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Obviously, thu tillipsoidal coordinate system (v,6,X) 
is orthogonal; and 0 and X are nothing but the so-called 
reduced latitude and longitude respectively. We then have 
an element of surface. 

(5) do - Gg 10 dX - Gg d:s 
an element of outer normal to o : 

(6) dn ■ Gjj dy 

where Gy , Gg and are called the first Gaussian funda- 

mental quantities defined bj. 

“ Et^'OSh^U- CO8^0]* 

°P *■ ‘ - E[cosh2y- cos*0]^ 

*^X * ^ ^ ^ ^ “ E coshycos0 “ Gj^COS 0 


It should be very easy to get their values on o : 
aO) G, = Gfl » c(i + e'2 sin0)^ + o(e'*) 

M P 

(11) ^ coshy s a 


Secondly, let us consider the Zagrebin's problem in the 
space (y,0,X) rather than (x,y,z). As a starting point we 
will take the potential T of the simple layer on o . With 
the action of the operator A the potential T , the outer 
limit value of the derivative of T in the direction of the 
normal on o , and the boundary condition which the potential 
must satisfy becomes respectively: 


(12) 

T ■ 

T" °0 ^x 

8 

ds 

PeCu 

(13) 

3 T 

-G- 

N 

-jr) Gg Gx ds - 2ir (^p G^p 

Peo 

(14) 

3T 

II 

- % 

Pco 


where s is a unit sphere, r , the distance from the para- 
meter point P to the integral variable point Q , and ^ 
means the density distribution function on the surface, Cu , 
the space outside the ellipsoid. 
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C(»nbinlng (12), (13) with (14), it follows that the Fredholm 
integral equation of the 2nd kind with respect to the function 
4 will be 


(15) - Ag + // 4K(p,Q)ds 

8 

where 

(16) K(P,Q) - g;> Gj 

is the kernel of the Integral equation. In our circumstance, 
it is also called the ellipsoidal kernel, because it is pro> 
duced fron the boundary surface of the ellipsoid. All quantities 
without the subscript P in (16) refer to the Integral variable 
point Q , but those outside the integral sign in (15) are 
always taken at the parameter point P on o . 

In order to solve the equation (15), it is necessary 
to represent the kernel (16) in more detail. The whole procedure 
is described as follows. 


It is clear that we have 


(17) 


and 

(18) 


^ ^ “ (x - Xp) t + (y - yp) 3 + (z - Zp) S 

r* » 4ppp sin*-|- + (p - pp)* 

p ■ (x*+y^+z*)^ •• E(8lnh*p + cos*B)^ ■ c(l + e' *cos*6) 


where ip is the angle between the vector radii ^ and Bp. 
Since the quantity p - pp is equal to 0(e'*): 


(19) p - Pp » i e' ^(cos*6 - cos*6p) 

the quantity (p-Pp)^ could be ignored, and (17) then becomes 

(20) r^ ■ 4ppp sin*-|- 


Substituting (18) into (20), we arrive at 

(21) r* « 4c*sln^>|- {1 + 4 e' *(cos*8 + cos*6p)} 
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(22) -I" ■ — T “T e'*(co8*0+ cos*^) ) 

2c8ln^ * ^ 

On the other hand, from (17) we have 

JiL. {rU p2 . p2 . 2(Z- Zp)Zp y 

3V^ ^ 8inh*i^ 


Using (21), (18), (4-3) and (4-4) it is not difficult to 
bring this derivative to the following form: 


(23) 4^ ■ rta^ , 114 1 { . (sine„ - sine)* ) 


3y 


Thus, we finally get 


4 sin*^ 


(24) 




** ** 4 ac sin -|- 


^ — r-{(l - i e‘*sin*8p - 


- i e'*(cos*0_ + cos*e) "*■ (siPB- - sin0)*} 

^ 4sin»-|- P 

We shall next derive the expression for the second tern 
in the bracket of (16). 


The Bruns' formula In our coordinate systen (y,0,X) should 
be the following form 


(25) 




W 1 1 

P P 



By means of the formulas of the radii of curvatui'e in the 
meridian and in the prime vertical 


(26) N • a(l e'*sin20)^ 

(27) U ~ ^(1 + -y e'*sin|> 
we immediately obtain 
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where 

(29) q - ^ 

Hence, we have 

{1 - e** sin* Bp + e’*sln*B q) 

Consequently, with (10), (11), (24) and (30) we find the 
expression of the kernal (16) as follows 

(31) K(P,Q) - K#(P,Q) + Ki(P,Q) 
where 

(32) Ko(P,Q) - ~r- 

4 sln-|' 

(33) Ki(P,Q) - e'* Ks{4-iisln*B„ + ^ sin*B-^ V -- - (sinB„ » 

o p 12 sin* ^ P 

“ slnB)* + ^q ) 

Moreover, with (22), (10) and (11), (12) may be written 
as 

(34) T - ^a // i(Kt ♦ k) d 8 

8 

where 

(35) k « e**Ko(i - i (co8*B - co8*6p)} 

We are going to determine the diaturblng potential T . 

For this purpose, we first solve the integral equation (15). 

For convenience, (15) is rewritten as the operator form 

(36) ^ ■ B4 "2^ 
where 

(37) Ut ■ Ag 




(30) 






**Yp 3Pp ac sin 
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. . . • « m * 


and ih9 operator B is defined by 

(38) B(. ) - /J K(* ) d 8 

Splitting (36) in the order of magnitude of e*’ without 
regard for o(e'^) we arrive at the following equations 

(39) 0(1): ■ Bo 08 ♦ ' 

(40) O(e'*):0i ■ Bo 0i 
where 

(40*) Gi - Bi 00 

and 

(41) B^; B^(«) - ^ // K^(*) ds 

is called the 1-th order operator with respect to e'^ . 

It should be easy to verify the following 3«»a: 

For any function X and G with the same order of 
e'^ , the operator equation ~ ~ 

(42) X - Bo X + 

always has a unique solution 

(43) 27tX - SG + G 

where S is the Stokes operator ; 

(44) S(0 - j-* // [S(0) - 11 (•) ds 

s 

Applying this iMBna to both equations (39) and (40) we have 

(45) 2w0o “ 5 Go ♦ Go 
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I 


and 

(46) 2iT(j)i » + Gi 

Consequently, we get the solution of the equation (36) 
as follows 

(47) ♦ = <|»o + ♦! 

Finally, we are going to detennine potential T . 

The expression (34) may be rewritten as 

(48) T = -I a[2TrBo<l> + ^ Gii] 

and also split into two parts in the order of e'* : 

(49) 0(1): To = a[2TTBo<|)o] 

(50) 0(e’^): Ti = a[2irBo4>i + ^ Gii] 

(51) Gii = J-a fj <|)okds 

s 

With (39), (45), (40) and (46) we get the solution of the 
Stokes' problem for the ellipsoid of rotation with a relative 
error of the flattening as follows 

(52) T = To + Tx 
with 

(53) To = S Go 

(54) Ti = S Gi + Gii 
or 

(55) 'To = ^ J/ Go [S(4;) - 1] d s , Go = Ag 

s 
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(56) Ti = -j- // GiESC'P) - 1] d s + -|-a // <l>o k d s , Gi “ ^ J/Ki4^,ds 

s s s 

Clearly, if the reduced coordinate system ($,X) is trans- 
formed into geocentric coordinate system (4»,X)» then our 
solution is reduced to the Koch's solution (Koch, 1968, p.22). 


4. Truncation Method of Rapp 


As mentioned above, the advantage of ellipsoidal kernel 
approaches consist in its closed iterative integral form. 

But in practical application, as we know, only numerical 
treatments can be carried out for these integrals. Moreover, 
they are extended over the global earth, so that the calcu- 
lation is fairly tedious and laborious even in order to obtain 
the geoid undulation of one point at the earth's surface. 

It is needless to say that it is much more difficult to eval- 
uate the geoid undulations of several points. Above all, 
a more important fact is that today people cannot have a 
uniform and dense distribution of gravity anomaly data on 
the surface of the earth. In this case, it is impossible 
to expect to arrive at good results through the procedures 
of global numerical integration. It follows that we have 
to find the other successful ways to reach our purpose. Even 
though the series methods are also complicated for practice, 
they bring light to the probability to use the information 
of potential coefficients to determine the geoid undulations. 
For containing the advantages of both series and ellipsoidal 
kernel approaches and overcoming their disadvantages a new 
ellipsoidal kernel procedure to compute the geoid undulation 
is developed by Lelgemann (1970). The principal part of 
his solution is exactly Stokes integral formula, and the 
correction term consists of a series of harmonic functions, 
unlike the general ellipsoidal kernel methods. Thus, we 
can use Lelgemann *s formula to ccMnpute the ellipsoidal cor- 
rection through the potential coefficient information. Later, 
a somewhat revised approach to the problem is suggested by 
Moritz (1980, p.320). 

It should be mentioned a well-known fact that the con- 
vergence of series of harmonic functions may not be very good. 
Theoretically, we must use a large number of terms to expect 
a rather good computed result. But in fact, the niimber of 
terms in the series of harmonic functions is only finite, 
in spite of either using gravity data measured at the earth's 
surface or combining these with the infoimation obtained 
by the aid of satellites. These coefficients cannot by any 
means represent the complicated situations of the gravity 
field of the earth. It follows that Lelgemann 's solution 
only takes into account the long wavelength part of correction 
tc the geoid undulation resulted from the earth's compression. 
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Aiming at this defect, Rapp (1981) developed another new ap- 
proach to calculate the ellipsoidal correction. He divided 
the solution of the problem into two parts: one is an inte- 

gral over a cap whose center is the computation point, which 
is evaluated by using gravity data measured at the earth's 
surface, and another is a series, calculated by means of 
the potential coefficient information. Thus, not only the 
effect of flattening on the long wavelength global undulation 
is considered, but also on the local undulation. This approach 
is preferable up to now from both theoretical and practical 
points of view. 

Rapp gives two kinds of formulas for the problem. We 
now write (Rapp, 1981): 

The first kind (using regular Stokes' function): 


(57) Ng = Ni + N2 + ANi + AN 
with 

R ,^02n 

(58) Ni = J J Ag SCip) sinij; dip da 

(59) N, = 2^ J Q„ 

(60) AN, = I (Q„-X„)Agi 

'' n=0 

(61) AN = e2(-^ - sin2(())N 


where 

(62) 

(63) 

(64) 


00 00 n 

Ag° “ I Ag° = I I 

n=2 n=2 m=0 

Ag>- I Ag; = -J- I (G 


n=0 


n=0 


H 1 

cos mX + sin mX )P 

n nm nm nm 

cos mX + sin mX )P ( sin(|) ) 

nm nm nm 


(sln(j) ) 


n 



n < 2 
n > 2 


For the meanings of the other symbols see (Rapp, 1981). 


As matter of fact, we can only calculate a finite sum 
of series; therefore we now give the truncation error of 



(65) 5Nj, « «Ni + 6 N 2 


where 



( 66 ) 

6 N 1 * ^ 1 

n=m+l 

(Q - X )Ag^ 
'^n n' *n 


00 


(67) 

6 N 2 = ^ 1 

n=m+l 

Q Ag» 
^n ®n 


According to (Hsu and Zhu, 1979), we have the variance 
of the truncation error of geoid undulation as follows: 

00 

( 68 ) (m,4'o) = (4)^ I (Ag«») 

6 n| n=m+l “ “ 


The second kind of computation (using the modified 


Stokes 

' function) is as follows: 


(69) 

N„ = Ni + N 2 + ANi + An 


where 

(70) 

R ,'^0,2'^ 

= 577/0 K ^etS( 4 <) - SirDo)] 

sinij^dtjfda 

(71) 

00 

^ i, 


(72) 

4N. =. ^ 1 (5„ - x„)4g‘ 

n«u 


(73) 

Qn “ Q„('i'o) + [Pjj_l(costi^o) 

-cos4>o Pjj(cost|;o )] , 


The truncation error of Nj, and its variance are ex- 
pressed respectively as: 


(74) 

6 Nj, =-■ 6 Ni + 6 N 2 


where 

(75) 

I , 

^ n=m+l 


(76) 

•SN 2 = 1 

^ n=»m+l 
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and 


/ R \2 


n-m+l 


5= o'(ig< 


5. Further Improvement of Convergence of Rapp's Method 


Rapp used the truncation theory to "remove" some local 
undulation terms from Stokes' series, the correction part 
in Lelgemann solution, into the integral over a cap surrounding 
the computation point, so that the convergence of the new 
series obtained is significantly accelerated, because there 
are not many local undulation terms in the series any more. 

We now suggest two procedures for further improving 
the convergence of Rapp's method. 

The first kind, using the optimum squares approximation 
to Stokes' function, is now described. 

We introduce the following approximation function 
S^(ip) to Stokes' function S(»|^) at interval [<Po ,Tr] : 

(79) » I ii^K_(m,ij;o) P„(cosij^) 

“ n=o ^ " 

where the coefficients Kjj(m,i|;o) satisfy the following 
condition: 

7T 

(80) sintj^dij.' = min 


From (Fang Tsiin, 1973, eq. (8.65)) we have: 

(81) Kn(m,ij^o) = Qn + I KCm,rl>o) 

s»0 ^ ns s 

< 1^0 

(82) X„c ~ J. P„(cosi|/) P_(costp) sirni'd^ 

II s 0 n ^ 


Let us analytically extend the function Sm(<|^) from 
region ['Po,it] to the whole region [0 ,tt] . Then, following 
(Rapp, 1981, eq. (40)), the solution of the problem can 
be expressed as 

“ "ify" 4 4 (Ag-e'AgM [S(4<) - S^,0P)] slnWda 
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“ e*Ag* )S_<ip)sinijjdt|)da 

^ >1 I 0 Q in 

+ j“ ^(Ag - e*Ag^rS(tp) - Sjjj(ip)]sim^di|^da + AN 

We write (Rapp, 1981, eq. (23)) 

(84) Ag = Ag® + e^Ag^ 


From (79) and (84) the second integral in (83) can be writ- 
ten as: 


(85) 


D 2 7T 

“ lu7 C 4 " e^AgMSjjj(\{») Sintfjdij^da- 



We are going to calculate the Integral involving e*Ag^ in 
the first term of (84). We can write: 


(86) sinij^dil^da 

= ^ - S^(ip)] sinij^dipda - 

" ^ /Jo sinil^dt^da 


We have: 


(87) 


sini^d,j;da 


Da 2 2 TT 

K L ^8 siru|;d^|^da « 

•» II T 0 0 


^ty C 4 Sinipdipda 

2 m 

‘ I 

n«o 


ly"nIo ^ „L Kn<®.4'o)Ag^ 

Re* ® 

n^O * 

+ T?” ^ ^^n “ ^n^^8n 

a-m+1 n n n 
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From (Hsu and Zhu, 1979) we arrive at; 

n«2 .■'f .21T , 

TrY ‘^'Po^o ) ]sint^)d\pda - 

= L, >''*i 

' n=m+i 

Hence, (86) becomes; 

R -»2 ,>JJo , 21 T 

(89) ^ AgMS(ii^) - S^(4^)]sin,|<di|;da 

p 2 ^ 

- !x„ - Q„ - + 

L tx„ - Q„ - Qi(m.t.)]Ag^ 

' n=m+l 

(90) Q;; = ^S(4-) - Sjjj(>p)] Pjj(cosiJ;) sinij;dijj 

Following (Hsu and Zhu, ibid) and (84) the third integral 
in (83) can be written as; 

(91) 6NJ = r - e^AgM [S(ip) - S (ij;)] sinij^dij^da 

ijj 0 ® *** 

7T 

“ I®**' ■ siiii|)di()da 

Therefore, the equation (83) finally becomes the fol- 
lowing form; 

(92) Ng = N{ + Ni + AN{ + 6Ng + AN 
where 

(93> N: = f'''Ag[S(i|;) - S„(iJ;)] sinij;d0a 

n ^ 

(85) N* = f- I K^(m,ij;o)Ag° 

n —2 

Pp 2 5 , 

(94) AN,: “ I Y (m,i|^o)AgJ 

n=0 " " 
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(95) 6N^ - 6Ni + 6N| 

(96) 6Nf - I Y'(m,,|^o)AgJ 

^ n*m+l 

(91) 6N1 = ^ I Q'(m.,j^o)Ag» 

n=m+l " " 

(97) Yjj(m,<»u> * K^(m,i|^o) + Q„ ’ Xn 

(98) Y;(m,4/o) “ Q^(m,ij;o) + Qn " X„ 

(61) AN = e*(-^ - -|-sln2(j,)N 


We have the variance of the truncation error: 

(99) ajj, .(m,i|;9)= (^)" I Q^^(m,4'o) aJ(Ag») 

E ' n=»ra+l 

From (Hsu and Zhu, ibid, sec. Ill) we find: 

(100) |Q^(m,^o)| < iQpl (and \^J) 

Comparing (99) with (68) and (77) respectively, we get 


(101) o^jj, (m,i|;o ) < o^jji(m,ii)o ) (and a^jjH(m,i|;o ) ) 

It means that the convergence of the series in case of the 
optimum approximation is better than that of the series 
in the truncation method of Rapp. 

The second kind (using the best squares approximation 
with a boundary restraint): 

In this circumstance, let the coefficients of approxima- 
tion function Sni('J^) the equation (79) be determined by 
the following minimum criterion (Hsu and Zhu, ibid, sec. V); 


(102) D = [S(ij;) - S^(i|;)]'siniJ^dif»do + 2Xi[S(ii^o) - S„(ii/o>] 

. 2X, 1^- -2^1- min 


- 16 - 


■ ffi 


in which, Xi and X 2 are undetermined multipliers. Using 
Lagrange's principle we arrive at approximation coefficients 

Kn •• 

(103) K - x''^ Q - (XiX"^U+ X2X"'W) 


(104) X - A“^B 


where 


K - [Ki, K 2 , . . . K^] 


X “ [Xi , X 2 ] 

Q 


[Qi I Q 2 » • • • Qjjjl 


u = 


u'^x'^u, u'^x’^w 
w'^x"'u, w'^x'^w 


iPo (cosif/o ) 

t 


2r+l 


P (costi^o ) 

I * 


2m+l 


Pm(cos4io ) 


B 


u'^x'^Q 

w'^x'^Q 


w 


2r+l 


r [P^(cosiJjo )-cos\poPj.(cosipo ) ] 


2m+l 


m.[P_(cosij^o )-cos\J;oP_(cosipo ) ] 

in in a 


Xo 0 » 

Xq 1 , . 

• 

• 

X 

0 

S 

1 

^•1 0 » 

Xi 1 , . 

• • ^im 

0 

e 

X 

1 

^ml’ • 

• • x^ 

mm _ 


>0 


X “/ P (cosip)P (cosip)sinipd({) 


We finally obtain the similar results to (92) and (99) 
but that the coefficients Kjj(m,\j^o) and Q'(m,ij^o) are now 
replaced with and Qi , where QJ5 may Be calculated 

by the following explicit approximate formula (Hsu and Zhu, 
ibid, eq. (50)): 
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(102) Q* - -.yl-— 8inij;oP_i(cos,|,,) Z(i)^,) + 0(— i— ) 

" n (n+1) n**/n 

(103) Z(\po) “ 2 + 9cosij;g + ^ | i -2S(i);o) 

/2 (1 - C08ll(o)^* 

It should be easy to see that |Qf| < |Q„| . It follows 
that the convergence of the series In this case is faster 
than that of the series of Rapp's method. 


6. Conclusions 


(a) All available procedures to approach Stokes' prob- 

lem for the ellipsoid as a reference surface could be clas- 
sified as two groups: the series solution, and the ellipsoidal 

kernel solution. 

(b) The ellipsoidal kernel solution represents a closed 
iterative form of Stokes' integral. It seems that this 

kind of solution is convenient for theoretical purposes. 

(c) The ellipsoidal kernel method developed by Lelgemann 
makes the computation of ellipsoidal correction to the geoid 
undulation by using potential coefficient information pos- 
sible; however, it can not be directly used for determining 
ellipsoidal correction to the geoid undulation through com- 
bining gravity anomaly information with potential coefficient 
information. 

(d) The method of Rapp, containing the advantage in 
Lelgemann method and removing its defect, is the best one 
of all available methods. It can be effectively used for 
determining highly accurate geoids from potential coefficient 
information and terrestrial gravity data. 

(e) As can be seen frcxn sec. 5, it is better to use 
the optimum squares approximation to Stokes' function in 
order to further improve the convergence of solution of the 
problem. 
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